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Abstract 

In this work we simulate the propagation of Ultra High Energy (UHE) 
protons in the magnetised intergalactic medium of Galaxy Clusters 
(GCs). Differently from previous works on the subject, we trace pro- 
ton trajectories in configurations of the Intra Cluster Magnetic Field 
(ICMF) which have been extracted from a constrained Magnetic-SPH 
simulation of the local universe. Such an approach allows us to take 
into account the effects of several features of the ICMFs, e.g. ir- 
regular geometrical structure and field fluctuations due to merger 
shocks, which cannot be investigated analitically or with usual nu- 
merical simulations. Furthermore, we are able to simulate a set of 
clusters which have properties quite similar to those of GCs observed 
in the nearby universe. We estimate the time that UHE protons 
take to get out of the clusters and found that in the energy range 
5 x 10 18 < E < 3 x 10 19 eV proton propagation takes place in the 
Bohm scattering diffusion regime passing smoothly to a small pitch 
angle diffusion regime at larger energies. We apply our results to 
estimate the secondary gamma and Hard X Ray (HXR) emissions 
produced by UHE protons in a rich GC. We show that the main emis- 
sion channel is due to the synchrotron HXR radiation of secondary 
electrons originated by proton photo-pair production scattering onto 
the CMB. This process may give rise to a detectable signal if a rela- 
tively powerful AGN, or a dead quasar, accelerating protons at UHEs 
is harboured by a rich GC in the local universe. 
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1 Introduction 



Several arguments suggest that Galaxy Clusters (GCs) may harbour sources 
of Ultra High Energy Cosmic Rays (UHECRs). If UHECRs are produced by 
sources located in the galaxies, e.g. compact stellar remnants or relativistic 
shocks produced by stellar collapse, these sources will be abundant in GCs 
which contains hundreds or thousands of galaxies. Active Galactic Nuclei 
(AGN), which are among the most promising source of UHECRs [1], are 
typically harboured by elliptical galaxies which are more frequently observed 
in GCs. Furthermore, GCs may act themselves as accelerators of cosmic 
rays up to ultra high energies due to the presence of large scale magnetic 
fields and shock waves produced during the cluster hierarchical accretion [2]. 

GCs are permeated by a magnetised plasma of electrons and baryons. 
The density of the Intra Cluster Medium (ICM) is determined by the ob- 
servation of the X-ray bremsstrahlung emission of the gravitationally heated 
electron gas. The inferred density is n gas ~ 1CT 3 cm~ 3 which amounts to 
an over-density of about one thousand with respect to the Inter Galactic 
Medium (IGM). The presence of Intra Cluster Magnetic Fields (ICMFs) is 
testified by the extended radio halos due to the synchrotron emission of rel- 
ativistic electrons in the ICM. The strength of the ICMFs can be estimated 
from the intensity of the observed radio emission either assuming the min- 
imum energy condition, giving < B >~ 0.1 — 1 /iG [3] (< B >~ 0.4 /iG 
for Coma [4]) or by an independent determination of the relativistic elec- 
trons density. For a few clusters (especially for Coma) this is made possible 
by the observation of a HXR emission that, if interpreted as the outcome 
of Inverse Compton Scattering (ICS) of the relativistic electrons onto the 
CMB photons, implies the ICMF strength to be in the range 0.2 — 0.4/iG 
[5, 6] However, the interpretation of the HXR radiation from GCs in terms 
of ICS is still controversial and other models have been proposed which may 
allow/require stronger ICMF (see e.g. [7, 8, 9]) . Indeed, Faraday Rotation 
Measurements (RMs) of polarised radio sources placed within the cluster, or 
in the background, provide significant evidence for the presence of stronger 
ICMFs in the range 1 — 10 fiG or higher [10, 11]. Furthermore, RMs provide 
valuable information about the spatial structure of the field which is patchy 
with coherence lengths in the range between 10 and 100 kpc. It is unclear 
if the discrepancy between the IGMF strength inferred from the radio halo 
intensity and that determined from RMs can be explained by the different 
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sensitivity of these two kind of measurements to the presence of magnetic 
substructures (for a discussion about this issue and a comprehensive review 
of IGMF observations see [12]). Interestingly, a recent analysis [13] of Bep- 
poSAX data on the HXR emission from Coma points to a ICMF strength 
in that cluster which is stronger (< B >~ 1 fiG) than the value previously 
claimed in [5] and it is almost compatible with that inferred from RMs. In the 
following we will assume that the actual strength of IGMFs is that inferred 
from RMs. 

The magnetised ICM can affect significantly the propagation of cosmic 
rays up to ultra high energies. The Larmor radius of protons in the ICMFs 
is smaller than the field coherence length up ultra high energies so that 
the residence time in a GC can be considerably increased respect to the 
straight propagation. The effect is even more pronounced for composite 
nuclei which have a smaller Larmor radius. Here we will focus on the case 
of protons. At low energies the propagation UHE protons is expected to 
take place in the spatial diffusion regime. Since diffusion arises due to the 
scattering of charged particles onto the field irregularities, the diffusion time 
is expected to depend significantly on the ICMF power spectrum. Although 
very interesting approaches have been recently proposed to determine the 
ICMF power spectrum from high resolution RMs [14, 15], so far this is a 
poorly known quantity. The common attitude is to assume a Kolmogorov 
power spectrum. Although this is the most natural choice in the case of fully 
developed, homogeneous turbulence of Navier-Stokes type (no-magnetic) we 
should take in mind that some, or any, of these conditions may not be fulfilled 
in the ICM. Indeed, the ICMF power spectrum may depend on the cluster 
accretion history. Furthermore, cluster merging is expected to compress and 
twist locally the ICMF giving rise to non-Gaussian MF fluctuations which 
are non accounted for in conventional CR diffusion simulations which use 
a synthetic MF with Gaussian fluctuations. The only way to account for 
the effects of all these ICMF features on the propagation of UHECRs is by 
means of numerical simulations of the MF large scale structure. Numerical 
simulations are also mandatory to simulate CR diffusion when the Larmor 
radius is of the same order of the MF coherence length. We will show that this 
is the case for protons with energy around 10 19 eV, i.e. the most interesting 
energy range for the extensive air shower (EAS) experiments. Furthermore 
the determination of the diffusion coefficients can be performed analytically 
only for MF in a regime of weak turbulence. We will show that this may not 
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be the case for ICMFs. 

The main aim of this paper is to estimate the diffusion time of UHE pro- 
tons in rich GCs as a function of their energy. Differently from previous works 
on the subject we approach this problem by means of numerical simulations 
to determine both the ICMF structure of nearby GCs and to trace proton 
propagation. We extract the ICMF in several GCs from a constrained simu- 
lation of the MF structure in the local universe [16] based on the Magnetic 
Smoothed Particle Hydrodynamics (MSPH) technique developed by Dolag 
et al. [17]. This simulation traces the passive evolution of the magnetic field 
in a spherical region of radius 115 Mpc about the Local Group starting from 
a primordial seed field. It was showed [16, 17] that the simulation repro- 
duces the observed RMs with good accuracy. The main advantage of this 
approach is that it provides a realistic simulation of the gas and of magnetic 
field spatial distribution in a number of nearby clusters. It was found that on 
average the ICMF traces the gas density and it is locally amplified in regions 
where accretion shocks fronts are present. The ICMF power spectrum has 
been determined down to the spatial resolution length which is ~ 10 kpc in 
the cluster center. Significant deviations were found in some GCs respect to 
a Kolmogorov spectrum. The effects produced by all these features of the 
ICMF onto the UHECR diffusion can be relevant and were not considered in 
the existing literature. Our approach allows to take them into account. 

As an application of our results we estimate the flux of HXR and gamma 
rays produced by the interaction of UHE protons with the gas and the radia- 
tion within GCs. We will mainly focus on the synchrotron emission of ultra- 
relativistic secondary electrons, produced by pp-scattering, which falls in the 
GeV region, and by proton photo-pair scattering onto the CMB falling in the 
HXR region of the electromagnetic spectrum. Whereas the former emission 
can hardly be detectable, the latter may be detectable if a relatively bright 
UHE proton source is harboured by a rich GC in the local supercluster. 

In Sec. 2 we give a brief description of the MSPH simulation and of the 
ICMF properties obtained from this simulation. In Sec. 3 we describe our 
ray-tracing code and the UHE proton diffusion properties at different ener- 
gies and radii. In Sec. 4 we show that the insertion of MHD turbulence 
at low length scales, which is unresolved in the MSPH simulation, does not 
affect significantly our results. In Sec. 5 we compare our approach with that 
followed by other authors. In Sec. 6 we estimate the secondary gamma-ray 
emission due to hadronic scattering onto the intra-cluster gas and the HXR 
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radiation due to the synchrotron emission of electrons and positrons gener- 
ated by proton photo-pair production scattering. Finally Sec. 7 contains 
our conclusions. 

2 MSPH simulations and ICMF properties 

The origin of ICMFs is still unknown. Although the observed high metallic- 
ity of the ICM suggests that it may have undergone a significant pollution 
driven by galactic winds , the strong intensity of ICMFs, their huge exten- 
sion (exceeding few Mpc's in some cases), and their large coherence length 
(10 100 kpc's) make hard to explain their origin uniquely by galactic ejec- 
tion. Another possibility is that ICMFs were generated starting from a seed 
which is subsequently amplified by the adiabatic compression and the shear 
flows of the gas driven by the hierarchical accretion of the clusters. Several 
mechanism have been proposed to explain the origin of the required seed field. 
For example, it could be ejected by AGN [19], or from a violent starburst 
activity at high redshift [20], it could be produced by some non-equilibrium 
process in the early universe [18], or be the result of a Biermann battery [21]. 
In principle, the battery has the advantage to be independent on unknown 
physics at high redshift since, in this case, the seed field is produced by 
thermoelectric currents powered by the cluster merger shocks. In practice, 
however, the battery can account only for too weak magnetic seeds [22]. A 
subsequent turbulent dynamo has to be invoked to increase the field strength 
enough to initiate a successful MHD amplification. Since this process cannot 
be simulated on the computer, the final intensity of the magnetic field is quite 
uncertain and has to be tuned to reproduce the RMs of GCs. Operatively, 
such an approach is, therefore, similar to assume a seed field of primordial 
origin. 

The simulation that we use here is based on the MSPH method developed 
in [23]. The code combines the merely gravitational interaction of the dom- 
inant dark-matter component with the MHD of the electron-baryon gas. In 
a previous work [17] it was showed that magnetic seed fields in the range of 
(0.2—5) x 10~ 9 G at redshift z* ~ 20 will be amplified due to the structure for- 
mation process and reproduce RM in clusters of galaxies. This corresponds 
to B = B(z 1f )(l + z*)~ 2 ~ (0.05-1) x 10~ n G at the present time in the un- 
clustered IGM. It was also demonstrated that the MF amplification process 
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completely erases any memory of the initial field configuration in high den- 
sity regions like galaxy clusters. Therefore, we can safely set the coherence 
length L c (z in ) of the initial seed field to be infinite in our simulation. 

In the simulation which we use here B = 10~ 12 Gauss. The initial 
conditions for the DM fluctuations were constructed from the IRAS 1.2-Jy 
galaxy survey by first smoothing the observed galaxy density field on a scale 
of 7 Mpc, evolving it linearly back in time, and then using it as a Gaussian 
constraint for an otherwise random realization of the ACDM cosmology. We 
extended the initial conditions of [24] by adding gas, together with an initial 
MF. Therefore the simulation at redshift zero represents the large scale mass 
distribution of our local universe as observed. Some of the prominent halos in 
the simulations can be identified with observed counterparts with mass and 
temperature within a factor of two or better. In [24] it is shown, that this 
simulation matches the observed local universe very well in many aspects. 
Therefore, such an approach allows to use a simulated set of GCs which have 
very similar properties to those of GCs observed nearby. 

The gravitational force resolution of the simulation is 10 kpc. This max- 
imal resolution is reached in the central region of GCs where the SPH test 
particles are more dense. The volume filled by the mixture of high resolution 
dark matter and gas particles is a sphere of radius ~ 115 Mpc (centered on 
the Milky Way) surrounded by a larger region of low resolution dark matter 
particles, to get long range tidal forces. The mass of the high resolution gas 
and dark matter particles thereby is 0.48 x 10 9 M e /h and 3.1 x 1O 9 M //i 
respectively. Therefore the most massive cluster in our simulation is resolved 
by nearly one million particles. 

For our analysis we used the four most massive clusters within our sim- 
ulation volume, which includes the halos identified as the Perseus and the 
Coma cluster. Table 1 summarizes the global properties of these clusters. 
We determined the gas density and magnetic field profiles for the clusters. 
On average, the magnetic field profiles follow the density profile in the outer 
part, whereas in the central region the magnetic field profiles flatten, which is 
consistent with previous findings. Depending on the cluster and its dynami- 
cal state, the slope of the magnetic profile in the outer parts scatters around 
the slope of the gas density, and also the size of the flatter core region varies. 
Fig. 1 shows, as an example, the profiles for the halo 2 and the halo 4. For 
the magnetic field we present the norm as well as the radial and transversal 
component. It is evident from these plots that the MF structure is basically 
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Table 1: Cluster number, virial radii, virial mass, emission weighted virial 
temperature and identification with real clusters (if possible). 

# R 200 [kpc/h] M 200 [M Q /h] T Lx [keV] name 



1 

2 
3 
4 



2421.76 
2257.75 
1958.90 
1885.23 



1.60e+15 
1.30e+15 
8.50e+14 
7.57e+14 



7.1 
5.6 
4.6 
8.3 



Perseus 



Coma 



isotropic. 




Figure 1: Density and magnetic field radial profiles of clusters 2 and 4. 
Different spatial components of the ICMF are represented. 

We also determined the MF power spectrum B 2 {k). As we will discuss 
in more details in the Sec. 3, this quantity is crucial to determine CR diffu- 
sion properties in the ICM. We computed B 2 {k) by evaluating the MF at 
the points of a cubic grid. This is done by applying the MSPH formalism 
to every mesh point of the grid considering all particles which overlap by 
their smoothing length, taken from the MSPH simulation. Afterwards we 
performed a FFT (Fast Fourier Transform) on this grid. 

In Fig. 2 we show the energy spectrum (B(k) 2 k 2 ) for clusters 2 and 4. The 
slope of the power spectrum for the four most massive haloes seems to be 
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steeper for the more massive ones, but when calculating the spectra for all 
clusters it looks more likely that it is determined by the dynamical state of 
the cluster. The spectral index 5 ranges between -1.5 and -3 for most of the 
clusters, with a small number of clusters found to have even more extreme 
values. On average we found 5 ~ —2 which corresponds to a spectrum slightly 
steeper than a Kolmogorov's for which (in 3D) 5 = —5/3. This is consistent 
with the results obtained in the first attempts to measure the ICMF power 
spectra from RM images of extended radio sources ([14, 15]). 




Figure 2: The energy power spectrum is represented for clusters 2 (left) and 
4 (right) in arbitrary units. The quoted values of 5 represent the best fit to 
the power law index. \5\ = 5/3 for a Kolmogorov spectrum. 

In order to quantify the turbulence strength we define the parameter 

(5B 2 ) _ ^J k k -dkeB\k) 
1 (B 2 ) jB 2 (x)d?x y } 

where k m i n is the wavenumber at which the energy spectrum gets its maximal 
value and the turbulent cascade sets-in. This is the scale at which most of the 
MF energy is concentrated. A; max is determined by the spatial resolution of the 
MSPH simulation. Although the physical wavenumber at which dissipation 
takes place is expected to be much larger, this artificial cut will not affect 
significantly r\ due to the rapid decreasing of the turbulent power at small 
scales. 
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From the Figs. 2 we can estimate /c m ; n ~ 50 ^ almost independently on 
the cluster mass. We have then a turbulence strength 77 ~ 0.6 over 1 Mpc 
box. In other words, our MSPH simulation predicts a quite high level of 
turbulence in rich GCs. 

3 Ray tracing and UHE proton distribution 
in galaxy clusters 

We simulate proton trajectories in the synthetic ICMF by solving equation 
of motion by means of a Runge-Kutta adaptive step size method. Our sim- 
ulation conserves energy with a very good accuracy: we get ^ ~ 10~ 5 in 
the worst case of 10 18 eV protons over a complete trajectory. As far as it 
concerns the simulation of proton trajectories, we disregard energy losses. 
We will verify the validity of this assumption a posteriori by comparing the 
residence time in clusters with the energy loss time due to the relevant col- 
lision processes. The ICMF is determined along the proton trajectory by 
performing a weighted sum (see [23] for details) of the magnetic fields of 
all smoothed gas particles in the MSPH simulation which overlap with the 
proton position. In same cases an extra turbulent MF component is artifi- 
cially added to the result of the simulation to model a possible unresolved, 
or physically unaccounted, small scale component of the IGMF (see Sec. 5). 

In Figs. (3,4) we represent typical proton trajectories in two and three 
dimensions. These plots clearly show that in the central regions of a rich 
GC proton propagation takes place in the spatial diffusion regime at energies 
as large as ~ 10 19 eV. Several momentum reversals are visible which we 
interpret to be due to magnetic mirroring onto local enhancements of the 
magnetic field strength. 

We start computing the mean delay time, i.e. the effective arrival time 
minus the straight propagation time to reach a sphere far out of the cluster 
(we take it at 35 Mpc radial distance), of 1000 protons injected at the cluster 
center with energies in the range 10 18 — 10 20 eV. In Fig. 5 we represent the 
delay time as a function of the proton energy for several GCs. The best fit 
function is well represented by two power laws with indexes ~ — 1 and ~ —2, 
the knee being located at an energy E cr ~ (2 — 3) x 10 19 eV. In Sec. 5 we will 
argue how this feature corresponds to the protons propagation passing from 
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5 x 10 18 eV proton in the central region of 



the spatial diffusion (SD) regime to the small pitch angle (SPA) regime at 
high energies. 

The delay time, that we defined in the above, is representative of the 
diffusive properties of the entire cluster. In order to get an insight on the 
radial dependence of this quantity, we computed the delay time from the 
center up to spheres of different radii. This is represented in the Figs. 6 at 
several energies for clusters 2 and 4. These plots show that no significant 
delay is produced at a radius larger than 1 Mpc . Interestingly, this feature 
is almost independent on the proton energy as it is apparent from the self- 
similarity of the curves. 

Our next step is to compute a set of time intervals Uj representing the 
mean time that a proton injected in the spherical shell Ri < R < Ri+i spend 
in the shell Rj < R < Rj+\ ■ Once we specified the injection spectrum and 
the radial distribution of sources, the knowledge of a sufficiently dense set of 
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Figure 4: The same of Fig. 3 in 3D and for protons with energy 10 19 eV. 

tij 's allows to determine the differential proton flux at different radii: 

dF(E,R i ) _ c A (R 3 j+ i-RD f(R , t m 

Here f(Ri) represents the density of sources and Q S {E) is the proton injection 
spectrum of a single source (for the sake of simplicity we assume all the 
sources to be identical). We took m = 2000 corresponding to as many 
identical shells each of them is 1 kpc thick. The normalization constant 
A depends on the source injection power. Eq. 2 amounts to a numerical 
solution of the diffusion equation. 

In the case of a point-like single source placed at the center of the cluster, 
f(Ri) = Sij, where 5ij is the Kronecker function. It is interesting to plot 
the histogram of t\j. For a better readability we represent it in Fig. 7 in the 
continuous limit. 

In Fig. 8 we represent the radial dependence of the proton energy density 
per logarithmic energy interval at several energies in the case of a single 
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Figure 5: The delay time (residence time within a sphere of radius 0.5 Mpc 
minus the straight propagation time) for 4 GCs. The quoted GC masses are 
in solar masse units. 



source, with luminosity L P (E > 10 18 eV) = 10 44 erg s -1 and spectral index 
a = 2 placed at the center of cluster 4 (Coma). 

Another possibility is that the UHECR source distribution traces that of 
galaxies. We can reasonably assume that the radial distribution of galaxies 
is well approximated by that of the gas. To be consistent we determine the 
gas distribution from the result of the MSPH simulation. We found that 




with r\ ~ lOkpc (spatial resolution limit), r2 — 250 kpc (core radius), r% ~ 
2 Mpc (virial radius) provides an excellent fit of the simulated gas distribution 
for GCs with mass M ~ 10 15 M Q . The mean UHECR luminosity of galaxies 
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Figure 6: The mean delay time accumulated from the cluster center to the 
radius R is represented as a function of R for protons of different energies 
(quoted in eV's) for clusters 2 (left) and 4 (right). 




Figure 7: The mean delay time accumulated in spherical shells of thickness 
1 kpc and mean radius R for clusters 2 (left) and 4 (right). 



can be estimated by assuming that the Milky Way provides a representative 
sample. The high energy cosmic ray luminosity of our Galaxy have been 
estimated by several authors, see e.g. [1] where it was found (we assume here 
that most of the UHECRs are protons Lf\E > 5 x 10 18 eV) ~ 5 x 10 36 erg/s. 
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100 1000 
R (kpc) 

Figure 8: The mean proton energy density in the cluster 4 (Coma) is plotted 
versus the radial distance from the center for several energies. 

Therefore, if the entire luminosity of a rich GC is due to ~ 1000 galaxies, 
its expected value would be in the range L p (E > 5 x 10 18 eV) ~ 10 39 — 
10 40 erg/s. We will show in Sec. 6 that this luminosity is too low to give rise 
to a detectable secondary emission from the ICM. 

A similar situation is expected if Gamma Ray Bursts (GRBs) are the 
main source of UHECRs. One of the main arguments in favour of this hy- 
pothesis [25, 26] is the coincidence of the predicted UHECR flux with that 
observed under the assumption that the energy release under the form of 
UHECRs by GRBs equals that of gamma rays. The rate of GRBs in a typ- 
ical galaxy is ~ 10~ 8 year -1 , hence it is Tq^ b ~ 10~ 5 year -1 in a rich 
GC. Our previous results imply that the diffusion time of UHE nuclei emitted 
in the central region of a rich cluster (R < 1 Mpc) exceed largely the time 
interval between GRBs in the cluster. As a consequence, clusters behaves as 
continuous sources of UHECRs even if elementary sources are GRBs. The 
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total UHECR luminosity of a galaxy due to GRBs is 

<b (E>5x «f eV) ~ (^j) = 10» ergs- . (4) 

This luminosity almost coincide with that previously estimated for UHECRs 
emitted by galaxies so that, even in this case, we do not expect an observable 
secondary emission from the ICM. 

4 Insertion of synthetic turbulence at small 
length scale 

As we discussed in the above, the maximal resolution of the MSPH simulation 
is about 2n/k m i n ~ 10 kpc. This means that the simulation does not account 
for magnetic fluctuations possibly present below this scale. Beside to the 
turbulent cascade, MHD turbulence at small scales may be due to a number 
of effects, e.g. galactic winds and galaxy motion. Since in the central region 
of GCs the Larmor radius of protons with energy below 10 19 eV is smaller 
than the spatial resolution of our simulation, the presence of strong magnetic 
fluctuations below this scale might spoil the validity of our previous results 
at low energies. 

In order to clarify this issue, we redetermined the delay time in the MSPH 
magnetic structure of a massive cluster to which we added turbulent MFs 
at small scale. The extra component of the MF is generated by summing 
over a large number of randomly distributed plane waves with spherically 
symmetric direction and with random polarizations and phases. We follow 
here the approach presented in [27] where it was showed that in the limit of 
an infinite number of wave modes, the turbulence is isotropic and spatially 
homogeneous. The modelled turbulent component is given by 

SB(x, y,z) = J2 A i k n) In exp(ik n z' n + ij3 n ) (5) 

where A(k n ) represents the amplitude of the wave mode with wave number 
k n , polarization £ n , and phase f3 n . Eq. (5) satisfies VB(i, y,z) = 0. The 
maximal value of k n that we used is 27r/0.5 kpc -1 . For the sake of simplicity 
we assume the extra component to have a Kolmogorov spectrum. Since the 
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MSPH simulation predict a stepper spectrum, we think that this assumption 
may only lead to overestimate the turbulent power at small scales hence the 
correction to the delay time. 

For generality, we assume the extra turbulent power to have the following 
radial dependence 

5B 2 (r) /r\- 

(6) 



Bmsph^o) Vr G 

Indeed, it is reasonable to assume that the level of turbulence will be higher 
in the center of GCs than outside. For example, since the formation of GCs 
the galaxies at the center make several revolutions around and may then 
contribute to the level of MHD turbulence. In this resonable choice 

for r is the cluster core radius r ~ 500 kpc. We adopt this choice to 
compute the delay time of protons versus their energy. This is plotted in 
Fig. 9 for different values of the parameter a. It is evident from this figure 
that the addition of turbulence at small scales doesn't increase significantly 
the diffusion time. This means that diffusion properties of protons in the 
energy range 10 18 — 10 20 eV are dominated by the MF fluctuations at larger 
scale, where most of the power is concentrated. Since these fluctuations are 
well accounted by the MSPH simulation we can reasonably trust the results 
presented in Sec. 3. 



5 Diffusion of UHE protons, comparison with 
previous results 



Diffusion of cosmic rays in disorganised magnetic fields takes place by scat- 
tering of the charged particles onto the magnetic irregularities. Particles 
interact with the field when their gyro-motion resonates with a Fourier com- 
ponent of wavelength equal to the Larmor radius. The diffusion coefficient, 

(Ax 2 ) 

defined by D(E) = , can be computed analytically in the regime of 

weak turbulence and for low values of the magnetic rigidity p = lixr^jLc <C 1. 
In this case 
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Figure 9: The same as Fig. 5 with an extra MHD turbulent power for several 
values of the parameters r and a defined in Eq.6. 



where B 2 (k)k 2 is the magnetic field energy power spectrum. In the case of 
a Kolmogorov spectrum, B 2 (k)k 2 oc k~ 5 ^ 3 , the energy dependence D(E) oc 
E~ l l 3 is obtained. A priori there are no reasons why Eq.(7) should hold for 
UHE protons diffusing in the ICMFs. In Sec. 2 we showed that MSPH predict 
strong turbulence (rj = 0.6) in the core of rich clusters (R < 500 kpc) with 
a power spectrum which may be significantly steeper than a Kolmogorov's. 
Furthermore, since the Larmor radius 

is comparable to the field coherence length L c , the magnetic rigidity of UHE 
protons in the cluster is quite close to unity. 

Numerical simulations of charged particle diffusion have been performed 
in the case of strong Kolmogorov turbulence [28]. These simulations con- 
firmed the validity of Eq.(7) for low values of the magnetic rigidity. However, 
significant deviations from the E~ x ^ 3 behaviour were found when the mag- 
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netic rigidity approaches unity. The best fit to the numerical results gives 
m units of Mpc 2 /Myr): 



( e y/ 3 ( b 



) 



1/3 



L 



) 



2/3 



D(E) 



0.004 



U0 20 evJ {^G 



Mpc 



E < Q.1E. 



cr 



D{E) 



D(E) 



0.02 



0.03 




(9) 



In these expressions E = E CT corresponds to the condition p — 1. Accord- 
ingly to these results, the energy dependence of the diffusion time, tdiff(-E) = 
r 2 /6D(E), changes rather smoothly from a quasi-rectilinear regime at E ~ 
E CI , with tdis(E) oc E~ 7 / 3 , to the, so called, Bohm diffusion regime with 
tdis(E) oc E^ 1 , for 0.1E cr < E < E CT . The analytically predicted behaviour, 
tam(E) oc E" 1 ^, was found only at low energies, E < 0.1E cr . 

However, the properties of the magnetic field configuration assumed in 
[28] are still quite different from those of ICMF given by our MSPH simu- 
lation. Indeed, the mean field intensity, which was assumed to be uniform 
in [28], has a quite pronounced radial profile in all the simulated clusters. 
Furthermore, MSPH simulations predict a turbulent power spectrum which 
is often quite different from a Kolmogorov's. We think that these simulation 
provide a much more realistic picture of the real ICMF configuration. 

First of all we checked that our numerical code reproduces the diffusion 
coefficients found in [28] in the energy range I0 18 — 10 20 eV for strong ho- 
mogeneous turbulence. Remarkably, we found that our results are in good 
agreement with those of [28] also in the case the field configuration is that 
determined by the MSPH simulation. The energy at which the delay time 
passes from a ~ E~ 2 dependence, consistently with small pitch angle dif- 
fusion, to E^ 1 (Bohm scattering diffusion) is ~ 3 x 10 19 eV. This energy 
corresponds to a Larmor radius of ~ 30kpc, for (B) ~ 1 fiG, which almost 
coincide with the scale L min = 2ir/ k rnin at which we found that most of ICMF 
energy is concentrated. 
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6 Hard X and gamma-ray emissions 



In this section we apply our previous results to estimate the flux of high 
energy photons produced by the interaction of UHE protons with the gas and 
the radiation in the magnetised core of a rich GCs. Such a radiation may 
provide a signature of the possible presence of powerful UHECR sources in 
galaxy clusters. For the sake of simplicity we will assume that a GC harbour 
a single point-like UHE proton source placed at its center. We will only 
consider the contribution of protons in the energy range 5 x 10 18 < E < 3 x 
10 19 eV with injection spectrum n p (E)dE = E~ J dE 1 . 

We start estimating the flux of gamma-rays due to hard hadronic inter- 
actions (pp-scattering) of UHE protons with the baryon gas in the ICM. The 
differential photon production rate Q^(E) from the decay of secondary 7r°'s 
is given by [2]: 

Q,(E) - Y, L ^ E > (JL) 1 a pp n gas ct diS (E) (10) 

where Y 1 is the fraction of proton energy transferred to the photons and 
o- pp ~ 3 x 10~ 26 cm 2 is the pp-scattering cross section (we neglect here a weak 
energy dependence of this quantity). In Sec. 3 we found that, in the energy 
range 5 x 10 18 < E < 5 x 10 19 eV, the diffusion path length of protons in 
the core of a cluster with mass ~ 8 x 10 14 M (like Coma) is 

Therefore the gamma-ray energy emission rate per logarithmic energy inter- 
val of the primary protons is 

40 ( ^7 \ f ^gas \ / Lp (E > -E^mm) \ / E \ ^ ^ 



E> Ql[ E) * 10- [ jjX ) [w ^) - 4erg - ^— ) erg s 



-'mm 



(12) 

where E min = 5 x 10 18 eV. Since the mean free path of E 7 ~ 10 19 eV photons, 
due to 7 + 7cmb e + + e~, is much smaller than the cluster core radius [29], 

1 By injection spectrum we mean here the spectrum of protons which leave the accel- 
eration region reaching the ICM. This may differ from the acceleration spectrum due to 
energy losses close to the source. 
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practically the entire energy of the photons will be converted into electro- 
magnetic showers inside the cluster. The shower, however, cannot escape the 
cluster. This is promptly understood by comparing the synchrotron energy 
loss length of the electron-positron component of the shower, 



with the cluster size and the IC scattering length over which electrons and 
positrons might energise CMB photons, which are much larger. Therefore, all 
the energy of the photons produced by 7r° decay is converted into synchrotron 
photons. For monochromatic electrons, and for a uniform magnetic field, the 
photon synchrotron spectrum would peak at the energy 



with width 5u ~ u syQ . Secondary electrons (and positrons) from pp-scattering 
however, are spread over a range of frequencies which is much larger than Sou 
so that the spectral shape of the synchrotron emission is mainly determined 
(we neglect here a possible dependence on the IGMF power spectrum which 
will be investigate elsewhere) by the electron energy spectrum [30] hence, in 
turn, by the primary proton spectrum in the cluster. 

It is worth noticing here, that the synchrotron emission in the GeV region 
would be overwhelmed by the secondary emission produced by low energy 
protons if the UHE protons spectrum extends to lower energies with a power 
7 > 2 . There may be cases, however, in which the proton emission is peaked 
at ultra high energies (see below about possible scenarios which may give rise 
to this situation). 

If we assume that the proton energy spectrum in the ICM peaks at 
-^max ~ 10 19 eV the synchrotron emission will peak at an energy of few 
tens of GeV's. Photons with this energy can travel over cosmological dis- 
tances without undergoing significant energy losses and might be detectable 
by the GLAST satellite [31] or the MAGIC Cherenkov telescope [32]. We 
find, however, that the expected photon flux 
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is much lower than the expected sensitivity of these instruments, unless the 
UHE proton source is extremely bright (L p (E > 5 x 10 18 eV) ^> 10 46 erg s _1 ) 
which, however, may contradict other constraints. 

A more promising signature of the presence of a powerful UHE proton 
source in a GC could be given by the synchrotron emission of secondary 
electrons and positrons produced by the process p + 7cmb e+ + e~ + P 
(proton photo-pair production). It is remarkable that, for protons in the 
energy range 5 x 10 18 <£<3x 10 19 eV, the secondary electrons produced 
by this process are practically monochromatic with energy ~ 830 TeV [1]. 
As a consequence, the spectrum of the synchrotron photons emitted in the 
cluster peaks in the hard X-ray (HXR) range, around ~ 10 keV (see Eq.14). 
Interestingly, a non-thermal emission from the Coma cluster, peaking at 
this energy, has been already detected. The observed energy flux is [5, 13] 
E\f(Ex) — 10~ n erg cm 2 s _1 . In order to estimate the primary proton 
flux required to explain this signal in terms of secondary e ± synchrotron 
emission we need to know the probability for a proton to undergo photo- 
pair production scattering inside the cluster. This probability is maximal at 
E ~ 10 19 eV, at which it takes the value 



where we used our Eq. (11) to find the proton diffusion time in the cluster 
and the interaction length l pa i r (E ~ 10 19 eV) ~ 1000 Mpc determined in 
[36, 37]. Therefore the UHE proton luminosity required to explain the Coma 
HXR emission entirely in terms of synchrotron emission of electrons produced 
by proton pair production scattering is given by 



This luminosity corresponds to L p (E > 1 GeV) > 10 46 erg s _1 for 7 > 2, 
and it increases rapidly for larger values of 7. Even for 7 ~ 2, the required 
proton luminosity would be too high to be compatible with the EGRET limit 
[33] on the secondary gamma-ray emission in the 100 Mev - 10 GeV range due 
to pp-scattering . It would also give rise to a too intense radio synchrotron 
emission frompp secondary electrons [7]. Furthermore, a proton source with 
luminosity L p (E > 5 x 10 18 eV) > 10 45 erg s' 1 may be at odd with the 




(16) 




(17) 
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results of UHECR experiments which, so far, do not show any evidence of a 
flux excess in the direction of Coma. 

A possible way out from some of these problems is to assume a flatter 
proton spectra in the cluster, i.e. 7 < 2. This may be the case if UHE protons 
are accelerated in the knots and hot spots of a powerful AGN's jet [34]. 
Another possibility is that the proton acceleration is not stochastic. Rather, 
it may be induced by the huge electric field generated by a supermassive 
rotating black- hole (a so called "dead quasar"). It was showed in [35] that 
proton energies as large as 5 x 10 19 eV can be reached in this case with a 
quite narrow spectral distribution. Although, a quite intense TeV gamma- 
ray direct emission is expected from such an object, this emission is predicted 
to be beamed so that its detection may be missed. In a such a case the X 
and gamma-ray diffuse emission produced by secondary particles in a GC 
harbouring the source may offer an independent signature. 

Even if the HXR excess from Coma is not due to an intense UHE protons 
emission in that GC, weaker proton sources may still give rise to a detectable 
HXR signal from closer GCs. Indeed the expected X-ray flux in the 10 — 
100 keV region due to synchrotron emission of secondary electrons is 

i2 /Ppair\ / d V (LAE > 5 x 10 18 eV)\ 2 

$ x ^2xl0" 12 — 35 1 ergs _1 cm~ 2 

V10- 2 / V20 Mpc/ v 10 43 ergs- 1 J S 

A proton luminosity L p (E > 5 x 10 18 eV) ~ 10 43 erg s 1 corresponds to 
L p (E > 1 GeV) ~ 10 44 erg s" 1 for 7 = 2.1. This is a standard luminos- 
ity for an AGN. The required luminosity can be even smaller if the proton 
spectrum is peaked at high energies. Therefore an AGN or a dead quasar 
accelerating protons at UHEs may give rise to a detectable signature under 
the form of synchrotron HXRs if they are harboured by a rich GC in the 
local supercluster. 

We note on passing, that a mechanism similar to that we just discussed 
here has been recently considered in [9]. The main difference with our mech- 
anism is that in [9] the synchrotron emitting electron-positron pairs were 
assumed to be produced by photons with energy E y = 700 TeV. No com- 
pelling process was suggested, however, to explain the origin of the required 
flux of photons with that particular energy. 

An independent signature of the presence of a UHE proton bright source 
may be provided by the secondary gamma ray radiation produced in the 
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IGM by the protons which escape the clusters. The relevant process at 
E > 10 19 eV is, again, proton pair-production onto the CMB. If IGMFs are 
weaker than 1CT 10 G, as it is suggested by the results of the large scale MSPH 
simulation performed in [16], the electromagnetic showers produced by the 
secondary electrons and photons can travel over a distance of hundreds Mpc's 
[38]. For a source at a distance of ~ 20 Mpc the electromagnetic shower at 
the observer position will be composed by photons, electrons and positrons 
with energies which can exceed several TeV's. Since the shower production 
probability is in this case P pa i r ~ 2 x 10~ 2 , the photon flux around 1 TeV 
can be roughly estimated to be 

A ^ ~ 1 TeV > = i^ p - = 2xir ' 2 (to^Ft) (20^) 2cmV1 • 

(W) 

Such a flux should not be missed by Cherenkov gamma-ray telescopes like 
MAGIC [32], HESS [39] and VERITAS [40]. 

Detailed computations of the expected spectra in the gamma as well as in 
the HXR region are beyond the aims of this work and they will be presented 
elsewhere. 

7 Conclusions 

In this paper we investigated the propagation of UHE protons in the magne- 
tised medium of Galaxy Clusters. By using a constrained MSPH simulation 
of the magnetic field structure in the local universe we were able to account, 
for the first time, for several features of the ICMF which we think to be 
present in actual nearby clusters. We showed that UHE propagation takes 
place in the spatial diffusion regime in the core of rich GCs (M ~ 10 15 M & ) 
up to the E CT ~ 3 x 10 19 eV. Below this energy, at least down to ~ 5 x 10 18 eV, 
the diffusion time scales with energy like E^ 1 (Bohm scattering diffusion). In 
spite of the different ICMF power spectrum found in different simulated clus- 
ters, we found that the behaviour of the proton residence time as a function 
of the energy does not changes significantly from cluster to cluster. Only the 
amount of the delay changes depending on the mean intensity of the ICMF. 
The path length increase due to the ICMF is too small to give rise to a sig- 
nificant modification of the proton energy spectrum due to energy losses in 
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the ICM, since ct dc ia y ~ 10 Mpc < l loss ~ 10 3 Mpc for E ~ 10 19 eV. Since 
the residence time of UHE protons is much larger than that expected for 
straight propagation, and for small pitch angle diffusion, the probability for 
them to undergo hadronic or photo-pair production scattering is consider- 
ably increased. In order to determine a possible signature of a bright proton 
source harboured by a GC, we applied our results to estimate the gamma 
ray and HXR secondary emission produced by UHE protons in the ICM. We 
showed that electromagnetic showers produced by pp-scattering and proton 
photo-pair production do not leave the cluster due to the intense synchrotron 
losses. Therefore, their energy is transferred to synchrotron photons of lower 
energy which can reach the observer without further losses. The synchrotron 
gamma-ray emission from the electrons produced by the decay of secondary 
pions is too weak to be detected. More promising is the synchrotron emission 
of secondary electrons and positrons produced by proton photo-pair scatter- 
ing which falls in the HXR range. We showed that the UHE proton emission 
of a relatively powerful AGN placed in a GC in the local supercluster may 
give rise to a detectable HXR emission. 
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